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Abstract 

The effect of mappings on the approximation, by Chebyshev collocation, of 
functions which exhibit localized regions of rapid variation is studied. A gen- 
eral strategy is introduced whereby mappings are adaptively constructed which 
map specified classes of rapidly varying functions into low order polynomials 
which can be accurately approximated by Chebyshev polynomial expansions. 
A particular family of mappings constructed in this way is tested on a variety 
of rapidly varying functions similar to those occurring in approximations. It is 
shown that the mapped function can be approximated much more accurately 
by Chebyshev polynomial approximations than in physical space or where map- 
pings constructed from other strategies are employed. 


Introduction. 

One of the major difficulties in the application of Chebyshev pseudo-spectral 
methods, or other spectral methods, to the solution of partial differential equations, 
is in the approximation of functions which exhibit localized regions of rapid variation. 
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The approximation of such functions by polynomial expansions generally results 
in global oscillations unless a large number of terms are used in the polynomial 
approximation. These oscillations often lead to instabilities and inaccuracies when 
spectral methods are applied to the solution of partial differential equations [8,11]. 
One possible remedy is to introduce a mapping to a new coordinate system so that 
the mapped function can be more accurately approximated by low order polynomial 
expansions. 

There is a great deal of computational evidence that appropriately chosen map- 
pings can significantly enhance the accuracy of pseudo-spectral applications. In [16] 
mappings were shown to enhance the accuracy of approximations to shock like (hy- 
perbolic tangent) functions. In [1,3, 5, 9] mappings were introduced dynamically via 
an adaptive procedure in which some error measure of the solution was minimized. 
Mappings have also been used to approximate boundary layer flows in Navier-Stokes 
calculations, for example [8]. 

There are two issues that must be addressed in the systematic application of 
mappings to enhance the accuracy of Chebyshev pseudo-spectral methods. These 
are 


• The construction of appropriate family of mappings 

• Criteria to choose a particular mapping from this family according to the 
behavior of the solution to be approximated. 

The second issue has been addressed by employing adaptive procedures in order to 
determine appropriate mappings [1,3, 5, 9]. In these procedures, a family of mappings 
is introduced depending on a small number of free parameters. In addition function- 
als are derived which are used to monitor the pseudo-spectral error. The appropriate 
mapping function is then chosen so that when the problem is transformed by the 
mapping the functional is minimized. 

The choice of an appropriate functional to monitor the pseudo-spectral error is 
an important component of the adaptive procedure. In [5,9] a Sobolev type norm of 
the solution was chosen as the functional monitoring the pseudo-spectral error. In [3] 
another functional, derived directly from the Chebyshev expansion of an arbitrary 
function, was used in an adaptive procedure and shown to be more effective than 
the Sobolev norm in measuring the errors of the pseudo-spectral approximation. 

Another equally important component of the adaptive procedure is the choice 
of the family of mappings to be utilized as coordinate transformations. Usually, 
the family of coordinate transformations to be incorporated in the solution method 
is prescribed in advance and justified by, at most, heuristic reasoning. There have 
been very few detailed studies of the effects of different families of mappings on the 
accuracy of pseudo-spectral approximations of various different types of functions. In 
[13] a family of mappings was derived which had the effect of giving a more uniform 
distribution to the collocation points, i.e. in an appropriate (but singular) limit 
the mapped Chebyshev collocation points were uniformly distributed. Numerical 
examples demonstrated that this mapping function gave enhanced accuracy over 
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the unmapped Chebyshev pseudo-spectral method for functions which did not have 
rapid variations. However the effect of the mapping on functions exhibiting rapid 
variation was not considered. A similar type of mapping was presented in [1]. In 
[3] a family of mappings was introduced which attempted to map certain interior 
regions to the boundaries. 

In this paper we derive a family of mappings which are well suited to the ap- 
proximation of functions with localized regions of rapid variations. The family of 
mappings depend on two parameters which are related to the size of the gradient 
over the narrow interval and the location of the region of rapid variation respectively. 
The effect of this mapping on the accuracy of a Chebyshev approximation is ana- 
lyzed for a variety of different functions and compared to the families of mappings 
in [3,13] and to the case where no mapping is employed. 

2.Pseudo-Spectral Approximation and Mappings 

We first describe the derivation of the Chebyshev pseudo-spectral approximation. 
This discussion is brief, a more detailed description may be found in [8,11]. Let f(x) 
be a defined in the interval 7, — 1 < x < 1. The Gauss-Lobatto points are defined as 

(1) x } = cos(^r), 0 <j<N. 

P n the polynomial of n th degree which interpolates / at the points Xj is the Cheby- 
shev interpolant, or pseudo-spectral, approximation to / and can be obtained by 
using the Gauss-Lobatto quadrature rule to evaluate the expansion coefficients. It 
is known that the maximum norm error in approximating / by the pseudo-spectral 
approximation differs by at most a factor O(logA) from the maximum norm error 
in approximating / by the minimax polynomial [14]. A similar result is true for the 
Galerkin polynomial approximation to /. We anticipate that the results presented 
here sire equally valid for Galerkin approximations. 

However, in general a polynomial expansion is not appropriate to approximate 
functions with large gradients, e.g. functions exhibiting spikes, localized oscillations, 
or near discontinuities. Furthermore, the accuracy of Chebyshev approximations 
tends to be sensitive to the location of regions with large gradients. There is com- 
putational evidence that Chebyshev approximations are more accurate for functions 
which vary rapidly near the boundaries as opposed to functions which vary away 
from the boundaries, for example [16]. It was shown in [16] that this is true both 
for Chebyshev methods, where the collocation points cluster near the boundaries, 
and also for collocation based on evenly spaced points. It appears from the results 
in [16] that this is a property of polynomial approximations and is not necessarily 
related to the clustering of the collocation points. 

We denote the family of mappings by, 

(2) x = q(s,a), 

where x represents the physical coordinate, — 1 < s < 1 is the transformed co- 
ordinate, and a denotes one or possibly more free parameters. The Chebyshev 
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interpolant can then be used to approximate the transformed function /(g(s,a)). 
The effect of the mapping can be regarded as transforming the function to be ap- 
proximated to /(g(s,a)) from /(x). If the mapping is properly chosen f(q(s,a)) 
will vary more gradually and so be more readily approximated by a polynomial. A 
different viewpoint is that in the original coordinate, x, we are changing the basis 
functions from polynomials to a new set of basis functions that can better represent 
functions that strongly vary in localized regions. 

The mappings can improve the accuracy of the pseudo-spectral approximation 
in three general ways which can be employed as strategies in designing families 
of mappings. The transformed function can be less rapidly varying in the new 
coordinate system so it is better approximated by a polynomial. Second the region 
of rapid variation can be mapped to the boundaries. Third the mapping can expand 
the region near the boundaries thus tending to provide a more uniform distribution 
of collocation points. This latter approach was proposed in [13]. The distinction 
between these three potential strategies is not sharp, for example a mapping in 
which a spike at x = 0 is transformed to a less rapidly varying function must of 
necessity map regions in the interior closer to the boundaries. However in this 
paper we will compare three families of mappings which were constructed based on 
these three potential strategies for several functions similar to those which occur in 
applications. Our results will demonstrate that in most instances mappings based 
on transforming the function into one which varies more gradually, tends to yield 
more accurate approximations. 

When the Chebyshev pseudo-spectral method is used to solve partial differential 
equations there are many sources of errors and it is difficult to differentiate the source 
of the errors. In this paper we concentrate on the interpolation problem and evaluate 
the effectiveness of different mappings in reducing the maximum norm error in the 
pseudo-spectral approximation. In such a program one could simply start with a 
function and find the mapping parameters which yield the smallest error in some 
norm. However in applications the exact solution is not known and some adaptive 
procedure is necessary in order to select appropriate mapping parameters. In this 
paper we will determine mapping parameters by minimizing a functional related to 
the spectral interpolation error. This functional, derived in [3], is 

(3) I 2 {g) = (^J^L 2 g) 2 * /w(s)ds S j , 

where 

w(s) = \/l — s 2 , L — w(s) — . 

ds 

It can be shown that (3) is equivalent to the standard Sobolev norm under the 
mapping x = cos(0). Furthermore, this functional gives an upper bound on the 
maximum norm of the error [3]. 

Our approach will be to compute (3) using the transformed function f(q(s,a)). 
We will do this by forming the Chebyshev interpolant to f(q(s , a)) using a fixed value 
of N and then, using the Gauss- Lobatto quadrature formula, evaluating (3). The 
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result will be a function only of the parameter vector a which we will denote by / 2 (a). 
We choose a so as to minimize / 2 . In practice this minimum is found by computing 
J 2 (a) over a range of values of a and then choosing that parameter for which a 
minimum occurs. In most cases we will take the global minimum over the range 
of values of a. In certain circumstances global minima occur for parameters which 
give a highly inaccurate approximation. This is due to inaccuracies in computing the 
quadrature in (3) when the transformed function f(q(s,a)) is poorly approximated 
for some parameters a. In this case we take local minima for which the Chebyshev 
interpolant is a good approximation to the given function /. 

This approach is a practical and effective method to determine mappings which 
enhance the accuracy of the pseudo-spectral approximation [3]. It has been success- 
fully utilized in the computation of problems in combustion [3,5,9] and has also been 
used for problems in wave propagation in [l] . Furthermore in all cases presented 
below the errors that were obtained from mappings for which /2(a) was minimized 
were very close to the minimum errors in the maximum norm over the range of 
parameters investigated. 

3.Families of Mappings 

We first introduce a family of mappings which have the property of transforming 
functions of a certain class to polynomials. Although the functions occuring in 
applications will not be exactly of the specified class, they would be expected to 
be sufficiently similar so that the Chebyshev approximation will require only low 
degree polynomials to give good accuracy. 

We begin with a specified class of functions of the form s = h(x , oi\ , 0L2) . We 
suppose for certain parameter values this function exhibits rapid variations. For 
example suppose that for large values of oq, this function exhibits shock-like behav- 
ior with the rapid variations occuring near the point x — 0L2 and nearly constant 
behavior away from this point. If the function h[x , oq, 0*2) is univalent then upon an 
appropriate rescaling this function can be assumed to map the interval — 1 < x < 1 
univalently onto itself. The inverse function x = c*i, <*2) then describes a 

family of mappings such that under this mapping the function h becomes linear and 
can therefore be approximated by a low order expansion of Chebyshev polynomials. 
In applications the rapidly varying solutions will not be exactly of the form of the 
given function h[x,a\, 0*2), however the image of these functions under the mapping 
given by /i -1 is likely to have a gradual variation and so can be approximated by 
a low order polynomial expansion provided the parameters a\ and a 2 are properly 
chosen. This will be justified by the examples presented below. 

The mapping x = h~ 1 (s, oc\, 012) can be expected to enhance the approximation 
not only rapidly varying functions similar to /i(x,ai,a 2 ) but also rapidly varying 
functions which behave like derivatives of h(x,OL\, a 2 ). In order to see why this is so 
let g = /i'(z,ai, a 2 ) where 1 denotes the derivative with respect to x. If /i(x,ai,a 2 ) 
behaves as a step function then g(x, oc \ , a 2 ) will behave as a spike centered at x = a 2 . 
Suppose that this point is located away from the boundaries. Using the chain rule 
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we obtain an expression for y(x,ai,a 2 ) under the mapping s = /i(x,ai,a 2 ), 

(4) g(h~ 1 (s,ai,a 2 ),ai,a 2 ) = (d/i -1 (s, a u a 2 )/ds)~ 1 . 

If the function /i(x,a 1 ,( 22 ) behaves as a step function changing rapidly near the 
interior point x = a 2 , then clearly dh~ 1 [s y ai, a 2 )/ds will change most rapidly near 
the boundaries in s. Therefore we expect the mapping to enhance the accuracy of 
the approximation for g. Similar arguments hold for the approximation of functions 
similar to higher derivatives of /i(x,ai,a2)* 

In this paper we consider the following function /i(x,ai,a 2 ) 

(5) s= s 0 + tan -1 (ai(x - a 2 ))/A. 

For large values of a\ this function is nearly discontinuous with a region of rapid 
variation occuring near x = c* 2 . The parameters so and A axe determined so that 

(5) maps the interval I onto itself. The values of so and A, determined in this way 
are 

So = — — , k = tan -1 (ai(l + a 2 ))/tan _1 (ai(l - a 2 )), 

K> -p J. 

A = tan -1 (ai(l - a 2 ))/( 1 - «o)- 
With these choices of A and so the inverse of (5), 

(6) x = a 2 + tan((s - s 0 )A)/ai, 

describes a two-parameter family of mappings of I into itself which is suitable for the 
accurate resolution of functions with localized regions of rapid variation. Here a\ is 
related to the degree of change of the function and 0:2 is related to the location of the 
region of rapid variation. In applications these parameters would would be obtained 
either from prior knowledge of the solution or from minimizing a functional such as 
(3) which measures the error in the approximation. We note that (6) is explicitly 
invertible. 

For small c*i, s is approximately equal to x and (6) is approximately the identity 
map. The use of the mapping (6) should not therefore be expected to degrade the 
accuracy in approximating slowly varying functions which can be well approximated 
by a Chebyshev approximation without any mapping. 

We point out that the an alternative to the inverse tangent function is the hy- 
perbolic tangent function with an appropriate scaling. Mappings constructed using 
the hyperbolic tangent are related to stationary shock-like solutions of the Burgers 
equation [6]. 

An alternative family of mappings has been proposed in [13]. This family of 
mappings is 

(7) x = arcsin(ais)/arcsin(ai) 

and 0 < c*i < 1. The effect of this mapping is to expand the boundary regions 
and compress the interior regions. As ai — ► 1 the image of the Gauss-Lobatto 
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points become uniformly spaced however the mapping becomes singular. As ct\ — ► 0 

(7) approaches the identity mapping. The analysis and examples presented in [13] 
demonstrate that the use of (7) with oq near 1 can significantly improve the accuracy 
of pseudo-spectral approximations for a range of functions. This was attributed 
to fact that the collocation points become nearly uniformly spaced as ai — ► 1. 
The major effect of this mapping is to obtain a more uniform distribution of the 
collocation points and as such the mapping depends on only one parameter (however 
a two-parameter extension was proposed in [13]). A related mapping is presented in 
[1] which permits a concentration of collocation points in the interior which exceeds 
that at the boundary, i.e. it goes beyond obtaining a uniform distribution of the 
collocation points. 

A third family of mappings that we will consider is 

(8) s = 4/n tan'^aitan^/^x' — 1)] + 1 
with 

x 1 = (a 2 ~ x)/(a 2 x — 1). 

The use of this mapping is described in [3,5]. The effect of the mapping is to map 
an arbitrary point, a 2 , to the origin and then expand one of the boundaries at the 
expense of the other by varying ai. This mapping is presented in [5] where it was 
used with a 2 = 0 to compute functions with rapid variations near the boundary and 
in [3] to compute functions with rapid variations in the interior. 


4.Mappings and Accuracy 

In this section we compare the effectiveness of the families of mappings described 
by (6), (7), and (8) in enhancing the accuracy of pseudo-spectral approximations. 
Our approach is to construct the pseudo-spectral polynomial for the transformed 
function and then to measure the maximum norm of the error. The maximum norm 
of the error is computed by comparing the approximating polynomial and the given 
function over a large grid of points. We also compute the discrete L 2 (unweighted) 
norm of the error in the coordinate system selected. We note that the L 2 norm of 
the error depends on the coordinate system and may not adequately measure the 
relevant errors in applications. 

We first consider approximations to the function 

(9) y = tanh (cr(x - x 0 )). 

with a = 40. In our computations we compute over a grid of values of ai and a 2 and 
select those parameter values for which the functional / 2 is minimized. For the case 
of (6) the global minimum always occurs when a 2 = xo, however, we indicate the 
sensitivity of the approximation to the value of a 2 . The approximating polynomial 
is computed using N = 121. The error is computed using the Gauss-Lobatto points 
with N = 351. Both the maximum norm and the discrete L 2 norm of the error are 
computed. The results are presented in Table 1. An entry of U in the first column 
denotes that no mapping was used. 


7 


Mapping 


a 2 

h 

Max Error 

L 2 Error 

6 

51.30 

0. 

2.40 

1.71e-14 

1.10e-14 

6 

6.06 

.2 

80.57 

8.37e-05 

2.18e-05 

7 

.99999 

NA 

132.54 

1.51e-03 

3.48e-04 

7 

.9 

NA 

181.24 

6.72e-03 

1.42e-03 

8 

9.92 

-.887 

189.12 

7.92e-03 

1.70e-03 

U 

NA 

NA 

227.14 

2.05e-02 

4.02e-03 


Table 1: y = tanh(<r(x — xq)) 


The data in this table illustrates the effectiveness of the mapping (6) in enhancing 
the approximation of the shock-like profile (9). A graphical illustration of this is 
shown in Figures la and lb. In Figure la we plot (9) against the original independent 
variable x. In Figure lb we plot the transformed function against the variable s using 
ai = 51.30 and a 2 = 0. The figures clearly illustrate the more gradual variation of 
the transformed function. 

We note that the mapping (7) allows for a more accurate approximation than 
if no mapping at all is employed and also permits a considerably more accurate 
approximation than the use of (8) . The best accuracy (and the minimum value of (3)) 
occurs as oq — ► 1 so that the collocation points are uniformly spaced but the mapping 
is singular. Values of a\ close to 1 do not affect our approximation program, however 
when mappings are applied to the solution of partial differential equations, Jacobians 
are introduced into the equations. The singularity in the mapping (7) which occurs 
at ai = 1 may therefore effect the solution of partial differential equations when 
c*i is sufficiently close to 1. In the solution of partial differential equations the 
degree of ‘near singularity’ that can be allowed in the mapping without affecting the 
numerical approximation will be problem dependent, depending on other parameters 
of the solution process not evaluated here. There is a degradation in accuracy as a\ 
is reduced from the singular value 1. We have presented results for ai = .9 in order 
to illustrate the behavior of the approximation as a\ is varied. 

Finally we note that (6) is robust in the sense that even if a 2 is not at the point 
of maximum variation it still yields a significantly improved approximation. This is 
important for several reasons. 

• The minimization might be inexact. 

• There may be several regions of rapid variation clustered together. 

• Different dependent variables might exhibit rapid spatial variation at different 
points. 

• In many two dimensional problems there is often a principal direction in which 
the spike (or rapid variation) occurs, however the location of this variation 
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Mapping 


a 2 

h 

Max Error 

L 2 Error 

6 

51.30 

0. 

2.40 

2.20e-ll 

1.70e-ll 

6 

6.06 

.2 

80.95 

1.96e-03 

6.68e-04 

7 

.99999 

NA 

119.22 

1.60e-02 

4.08e-03 

8 

9.92 

-.887 

137.60 

5.11e-02 

1.12e-02 

U 

NA 

NA 

130.77 

8.68e-02 

1.92e-02 


Table 2: y = tanh(a(x — xo)), using N = 81 

may depend on some transverse coordinate. It may be more efficient to em- 
ploy a one-dimensional mapping to resolve such a solution rather than a two- 
dimensional mapping which depends on the transverse coordinate. 

The behavior described in the last point above is typical of partial differential 
applications, for example in wave propagation if a pulse propagates through a dis- 
persive medium, in fluid dynamics if curved shocks form and in combustion where 
cells can form along flame fronts. Thus we consider the robustness of the mapping 
(6) with respect to the location parameter a 2 to be an important point. 

We next illustrate the convergence properties for the mapped function by con- 
sidering the same results using N = 81. The results are presented in Table 2. In 
this case we have used the same mapping parameters as were obtained using 121 
collocation points. The results show that there is a much more rapid convergence 
when the mapping (6) is used. This is probably due to infinite order convergence of 
the spectral method which in practice means a very high order of convergence when 
the function varies gradually. We note that there are definite numerical errors in the 
approximation of the function / 2 when the mapping (6) is not used. In practice this 
could lead to poor parameter values being produced by the minimization procedure. 

The results described for (9) are not necessarily surprising as the family of map- 
pings (6) is constructed so that the image of a near step function is linear. We next 
consider the function 

(10) y = exp(a 2 (x - x 0 ) 2 /2). 

This function also exhibits a rapid variation near x = xo however the behavior is 
that of an isolated spike rather than that of a near step as in (9). We consider the 
values g = 50 and xo = 0. The results are presented in Table 3. 

We can infer from the results presented in Table 3 that the mapping (6) is 
effective in enhancing the resolution of spike-like functions. We also illustrate this 
graphically in Figures 2a and 2b where we plot the function (10) in the original 
coordinate, x, (Figure 2a) and the transformed function in the mapped coordinate, 
s, (Figure 2b). We note that the mapped function still resembles a Gaussian, however 
there is a much more gradual variation. 

The results also illustrate that the use of the mapping (6) provides good accuracy 
as the parameter a 2 is moved away from the location of the spike. The results 
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Mapping 

1 

a 2 

h 

Max Error 

L 2 Error 

6 

37.48 

0. 

18.26 

1.42e-09 

8.30e-10 

6 

10.86 

.1 

106.67 

3.37e-06 

1.48e-06 

6 

5.744 

.2 

245.51 

1.59e-03 

5.51e-04 

6 

3.506 

.3 

386.86 

2.20e-02 

6.08e-03 

7 

.99999 

NA 

285.20 

1.20e-02 

6.06e-03 

7 

.9 

NA 

557.02 

4.90e-02 

1.96e-02 

8 

73.92 

-.983 

547.04 

4.57e-02 

1.85e-02 

U 

NA 

NA 

699.84 

1.14e-01 

3.95e-02 


Table 3: y = exp(<r 2 (x — xq) 2 /2), xq = 0 


suggest that the mapping (6) can allow high resolution of multiple spikes even when 
the parameter c*2 deviates significantly from the location of each individual spike. 
Of course if the spikes are widely separated or if ol<i is far from the center of the 
spike the effectiveness of (6) degrades. We have found that the accuracy obtained 
from using (6) does not strongly depend on the value of xq provided a 2 is chosen 
to be near xq. Similar results hold for functions exhibiting a step function type of 
behavior as in (9). We also note that the mapping (7) provides considerably more 
accuracy than (8) indicating that improvements over the original Chebyshev method 
when the boundary clustering of the collocation points is reduced. 

The results in Table 3 indicate that the functional (3) has a rather sharp mini- 
mum when 0 L 2 is at the location of the spike. In our calculations we have found the 
value of c*2 obtained from minimizing (3) to be an excellent indicator of the location 
of the spike (or most rapid variation). This suggests that this procedure could be 
used as a shock locator for non-oscillatory spectral methods and filtering methods 
[7,10], at least when only one shock is present or multiple shocks are closely spaced. 

The effect of the location of the spike does considerably effect the behavior 
obtained from using (7) and (8). We illustrate this in our next example where we 
consider a spike near the boundary by setting xo = .9. The results are presented in 
Table 4. 

The mapping (6) still gives very good accuracy even when <*2 is not close to 
xo- The best approximation using (7) now occurs when ol\ — ► 0 where the mapping 
approaches the identity. In Table 4 we indicate the errors found for a\ = .0001, 
the limits of our search region together with the errors for larger values of oc\ which 
lead to approximately uniformly spaced points. For functions which have significant 
variations near the boundaries having the collocation points more uniformly spaced 
leads to a degradation in accuracy. We find that the behavior of the mapping (7) 
depends very abruptly on the location of the spike. For values of xo just below .8 we 
find an abrupt transition where the value of ai yielding maximum accuracy switches 
from a 1 near 1 to near 0. In this narrow overlap region the results obtained from 
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Mapping 

OL 1 

OL 2 

h 

Max Error 

Li Error 

6 

35.25 

.9 

15.05 

3.29e-10 

2.21e-10 

6 

1.89 

.5 

179.55 

1.29e-04 

5.10e-05 

7 

.0001 

NA 

227.22 

1.01e-03 

4.25e-04 

7 

.9 

NA 

288.03 

2.85e-03 

6.43e-03 

7 

.99999 

NA 

399.82 

1.54e-02 

6.43e-03 

8 

..0896 

0 . 

54.97 

8.31e-10 

4.03e-10 

U 

NA 

NA 

227.22 

1.01e-03 

4.25e-04 


Table 4: y = exp(cr 2 (x — xo) 2 /2),xq = .9 


Mapping 

<*i 

a 2 

h 

Max Error 

1-2 Error 

6 

5.136 

.0738 

386.52 

2.75e-03 

1.37e-03 

7 

.99999 

NA 

567.2 

1.56e-02 

8.17e-03 

7 

.90000 

NA 

787.1 

7.81e-02 

2.60e-02 

U 

NA 

NA 

1031.0 

1.70e-01 

5.54e-02 


Table 5: Sum of two Gaussians, xo = .30, x\ = —.15 
(7) are insensitive to a\. 

The results using (8) were found to be very insensitive to c *2 and are presented 
only for the mapping with oli = 0 as it used in [5] for functions with rapid variation 
near the boundary. The results indicate that this mapping is effective in approximat- 
ing boundary layer type solutions however it does not seem to offer any advantages 
over (6) and is not considered further. 

We next consider the following function 

(11) y = exp(<r 2 (x - x 0 )) 2 /2) + exp(o- 2 (x - xi)) 2 /2) 

The objective will be to determine the behavior of the various mappings when there 
is more than one region of rapid variation. The separation of the two regions of 
rapid variation can be controlled by varying xo and x\. We first consider the case 
xo = .30, xi = —.15,0- = 50. The results are presented in Table 5. 

It can be seen from the table that the family of mappings described by (6) yields 
a more accurate Chebyshev approximation even in the case of relatively separated 
pulses. We note that the minimum error of the mapping (7) in this case again 
appears to occur when ol\ — ► 1 so that the collocation points are nearly uniformly 
spaced but the mapping is almost singular. We have also exhibited the degradation 
in accuracy obtained by reducing ol\ from 1 by specifying an upper limit on ol\ of 
0.9. Even with this degradation the results obtained from using (7) are considerably 
better than for the case when no mapping is employed. 
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Mapping 

oti 

C* 2 

h 

Max Error 

Z >2 Error 

6 

2.704 

.0780 

625.4 

4.73e-02 

1.43e-02 

7 

.99999 

NA 

567.1 

1.45e-02 

7.96e-03 

7 

.90000 

NA 

758.8 

5.76e-02 

2.35e-02 

U 

NA 

NA 

1018.04 

1.95e-01 

4.84e-02 


Table 6: Sum of two Gaussians, xq = .50, x\ = —.25 


Mapping 


a 2 

h 

Max Error 

L<i Error 

6 

3.172 

.361 

375.3 

3.86e-03 

1.57e-03 

7 

.99999 

NA 

566.6 

1.43e-02 

7.71e-03 

7 

.90000 

NA 

629.3 

8.01e-02 

1.99e-02 

U 

NA 

NA 

742.91 

1.44e-01 

3.78e-02 


Table 7: Sum of two Gaussians, xq = .90, x\ — .2 


We next consider more widely separated pulses by setting xo = .5,xi = —.25,0- = 
40. The results are presented in Table 6. This is a relatively severe case as neither of 
the pulses are located close to the boundary. We observe that for widely separated 
pulses located away from the boundaries the use of (6) leads to larger errors than 
(7), provided values of a\ near 1 are used. The results obtained from (7) are not 
strongly sensitive to xo. This is not surprising as the effect of the mapping is to 
make the collocation points more uniform (in x). 

We finally consider the case where one of the pulses is located near the boundary. 
In this case we consider xo = .9,xi = .2, a = 50. The results are presented in Table 
7. We note that the best approximation occurs when c *2 is close to the location of 
the interior spike. 

In the case of multiple regions of rapid variation the optimum values of a\ 
and c*2 for the mapping (6) are not suggested by the locations of the regions of 
rapid variation. These parameters must be determined by a numerical minimization 
procedure. 

In summary our results indicate that the strategy of transforming the given 
function to one that is more readily approximated by a polynomial, which led to 
the derivation of the family of mappings (6) appears to be more effective than the 
strategy of mapping points of rapid variation to the boundaries. In particular the 
use of (6) appears to be superior in all respects to (8). The parameter a<i in (6) 
is an excellent indicator of the location of maximum variation, but the mapping is 
very robust if c *2 is varied away from this location. 

When the rapid variation occurs away from the boundary the use of (7) with ol\ 
near 1 so that the collocation points are nearly uniformly spaced is preferable to the 
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Mapping 

1 

2 

h 

Max Error 

Li Error 

6 

6.54 

0 . 

3.61 

6.00e-9 

4.23e-9 

7 

.99624 

NA 

7.99 

3.66e-05 

3.06e-04 

U 

NA 

NA 

14.59 

7.80e-04 

5.36e-04 


Table 8: y = 0.05/(x 2 + 0.05) 

unmapped Chebyshev method. In these cases the mapping (7) is most accurate as 
ai — ► 1 and degrades when oq is reduced from the singular value. It is therefore 
preferable to use this mapping as close to the singular point as possible. When there 
is a single region of rapid variation near the boundaries the unmapped Chebyshev 
approximation is more accurate than the use of (7). 

The use of (6) is superior to (7) for functions with rapid variation except when 
there are multiple regions of rapid variations which are widely separated and not 
located near the boundaries. In this case the best strategy seems to be to use (7) to 
obtain nearly equally spaced collocation points. Although this analysis is presented 
for pulse like functions similar effects hold for other types of functions with regions of 
rapid variation, for example functions which exhibit a step function type of behavior 
or localized highly oscillatory behavior. 

In this paper we have concentrated on the role of mappings in enhancing the 
accuracy of Chebyshev pseudo-spectral approximations. The transformation (7) 
was originally introduced to decrease the spectral radius of the spectral differenti- 
ation operator and thereby increase the allowable timestep in the solution of time 
dependent partial differential equations. 

Although the major emphasis in this paper is on the approximation of rapidly 
varying functions we have found that the mapping (6) is also effective in approxi- 
mating functions which vary more gradually. We illustrate this with the example 

(12) y = 0.05/(x 2 + 0.05), 

which was also considered in [13]. In this case we use 32 collocation points. The 
results of the approximation are shown in Table 8. 

We next illustrate the effectiveness of this mapping in approximating solutions 
of partial differential equations. We consider the Burgers equation 

(13) u t + uu x = vu xz , 

where — 1 < x < 1. The initial condition and the boundary conditions are chosen as 

u(x,0) = — sin^rx, u(±l,t) = 0. 

The viscosity coefficient v = .01/tt. This problem has been used as test case for a 
variety of spectral, pseudo-spectral and finite difference methods [2,9]. As t increases 
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m 

'Kt 

ai 

151.9896 

1.60407 

21.5 

152.3214 

1.60404 

37.0 


Table 9: Burgers equation with N=61 

the solution develops a very steep gradient at the x = 0. The maximum gradient at 
x = 0 occurs when nt is close to 1.6. As t increases beyond this point the solution 
relaxes to zero. We have solved this problem using an adaptive pseudo-spectral 
method together with the mapping (6). The solution is symmetric about x = 0 and 
the minimum of (3) always occurs at cl^ — 0. 

The solution to this problem can be computed analytically. One measure of 
accuracy, which has been used in [2,9] is the quantity 

m = max(u z (0, t)). 

The analytic result quoted in [2] is m = 1.5200516 occuring at nt = 1.6037. We 
have found analogously to what was found in [9] that the coordinate system which 
yields the most accurate calculation of m is not perfectly predicted by the adaptive 
procedure. Minimizing the functional (3) appears to concentrate too much resolution 
near x = 0 at the expense of resolution away from x = 0. A similar conclusion 
was arrived at in [9] using an adaptive procedure based on minimizing a different 
functional. When the search region is artificially constricted so that the parameter 
c*i does not get too large, effectively reducing the resolution near x = 0. we find 
excellent agreement with the analytic value for m. We indicate the computed values 
of m in Table 9. The number of collocation points, A, was 61. The first entry in 
the table is the result when the value of a\ was restricted while the second entry is 
for the value of aq selected by the adaptive procedure. 

Using the mapping (6) a high degree of accuracy can be attained as the number 
of collocation points is reduced. We indicate this in table (10) where we list the 
values of m obtained for different numbers of collocation points. In this table we 
have not restricted the search for aq but rather taken the value of oq predicted 
by minimizing (3). In most calculations the exact solution is not available and the 
minimization of some functional is at present the only practical method to predict 
the appropriate coordinate system. 


4.Multiple Domains 

The results presented above indicate that the mapping (6) can significantly im- 
prove the Chebyshev approximation to rapidly varying function when there is only 
one region of rapid variation or when there are several closely spaced regions of rapid 
variation. When there are widely spaced regions of rapid variation away from the 
boundaries the effectiveness of (6) degrades and the best strategy within the con- 
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m 

nt 

OL\ 

N 

151.7508 

1.61622 

37.0 

33 

152.3011 

1.60446 

37.0 

49 

152.3214 

1.60404 

37.0 

61 


Table 10: Burgers equation varying N 

text of a single domain is to use (7) with a near 1 so as to obtain a nearly uniform 
distribution of collocation points. However the accuracy obtained from using (7) is 
generally less than the accuracy that could be obtained from use of the mapping (6) 
if the regions of rapid variation were not widely separated. 

One possible approach to this problem is to extend the mapping (6) to allow for 
the resolution of multiple regions of rapid variation. Such an extension would re- 
quire additional parameters in the mapping and therefore the resulting minimization 
problem to determine appropriate parameters would be considerably more expensive 
than for the case when the mapping depends on only two parameters. 

An alternative approach is to introduce two or more subdomains and employ 
mappings within each subdomain. The use of multiple domains can lead to an 
improvement in accuracy by (a) resolving small scale structures in the problem by 
introducing domains corresponding to the length scales appropriate to the problem, 
(b) choosing the interface so that the small scale structures or rapid variations occur 
near the boundary, and (c) isolating different regions of rapid variations within each 
subdomain and then employing mappings such as (6) within each subdomain in 
order to improve the accuracy of the unmapped Chebyshev approximation. 

Multiple domain procedures have other advantages in addition to possible im- 
provements in accuracy, for example they lend themselves to parallel computation 
and can lead to smaller and better conditioned matrices, see for example [8,15,12]. In 
this paper however, we will only consider the effect of multiple domains on the accu- 
racy of approximation of rapidly varying functions. In particular we will demonstrate 
that a strategy such as that described by point (c) above, involving the interaction 
of mappings and domain decomposition can lead to significant improvements in the 
accuracy of Chebyshev interpolation. 

In the context of the solution of partial differential equations, the use of multiple 
domains introduces addition al parameters which have to be determined, in particu- 
lar the number of domains, the location of the interface points, and the number of 
collocation points within each subdomain. If the multiple domains are chosen to cor- 
respond to very small scale structures or localized rapid variations, then the accuracy 
of the resulting approximation may be may be quite sensitive to the location and 
extent of the interface points. These parameters are generally chosen non- adaptively 
or with prior knowledge of the behavior of the function to be approximated. 

When mappings are employed additional parameters must be determined. For 
example, if there are two domains and the mapping (6) is employed in each sub- 
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domain then there are 7 parameters which have to be determined (the parameters 
of the mapping, the location of the interface and the number of collocation points 
within each subdomain). In this paper we will constrain the number of points so 
that each domain has an equal number of points. This strategy is preferable for 
parallel computation, although it may not necessarily be preferable for determining 
the most accurate approximation with a minimum number of computational degrees 
of freedom. We will also constrain the number of domains to be 2. We now require 
the determination of 5 parameters, together with some procedure to determine the 
location of the interface, e.g. some functional which measures the error as a func- 
tion of the location of the interface. The determination of these parameters could be 
significantly simplified if the location of the domain could be determined indepen- 
dently of the mapping parameters within each subdomain, i.e. a strategy whereby 
the location of the domain is determined from the current solution, and then new 
mapping parameters are determined within each subdomain by the minimization 
procedure described above. We will demonstrate that a strategy of employing mul- 
tiple domains to isolate regions of rapid variation, together with the mapping (6) to 
resolve localized regions of rapid variation within each subdomain can lead to such 
a decoupling of the problem of determining these parameters. 

In order to accomplish this, it is first necessary to understand the role of the 
location of the interface on the accuracy of the approximation and to develop func- 
tionals which can monitor this. We will therefore first consider the use of multiple 
domains without any mapping within each subdomain and then demonstrate that 
once an appropriate interface is found the errors in the approximation can be fur- 
ther reduced by employing mappings such as (6) within each subdomain. We will 
consider the approximation of (11) with a = 40 and with 41 points within each 
subdomain. We will further restrict ourselves to the case where xo = .5,xi = —.25 
although similar results have been obtained for other values of these parameters. 
The accuracy of single domain Chebyshev approximations to this function under 
various mappings is described above. We use 41 points within each subdomain and 
compute the largest of the maximum norm errors. 

The location of the interface, Q, will be obtained adaptively by determining Q 
so that the maximum of (3) in each subdomain is minimized. Since we use an equal 
number of points in each subdomain, this corresponds to minimizing an estimate 
for the largest of the maximum norm errors. In Table 11 we show that this strategy 
leads to an improvement in accuracy over the single domain result and in the case 
where the interface is chosen in an ad hoc manner to be equidistant from both 
spikes, i. e. Q = (xo + xi)/2. We have then employed the mapping (6) within each 
subdomain to determine parameters a},a 2 > a i> a 2 > determined to minimize (3) in 
each subdomain els in Section 4 above ( we employ a simple linear mapping to map 
each subdomain to the interval I). 

It is apparent from the Table 11 that the overall accuracy of the approximation 
is very sensitive to the location of the interface when no mapping is employed. In 
this case the strategy of minimizing the maximum, over the two domains, of (3) 
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H 

If 

Q 

Mapping 

a \ 



— 7? — 
a 2 

Max Error 

1018.04 

single domain 

u 

NA 

NA 

NA 

NA 

1.95e-01 

625.4 

single domain 

(6) 

2.7404 

.0780 

NA 

NA 

4.73e-02 

567.1 

single domain 

(7) 

.99999 

NA 

NA 

NA 

7.96e-03 

215.8 

159.1 

.125 

U 

NA 

NA 

NA 

NA 

1.08e-01 

17.52 

17.24 

.125 

(6) 

20.60 

.330 

15.70 

-.1417 

5.17e-06 

183.0 

181.8 

.05 

U 

NA 

NA 

NA 

NA 

6.44e-02 

17.25 

17.56 

.05 

(6) 

47.75 

.425 

17.47 

-.0473 

4.59e-06 

150.27 

131.05 

.05 

(7) 

.99999 

NA 

.99999 

NA 

3.29e-02 


Table 11: Effect of mappings on multi-domain approximation 


leads to a significant reduction of the error. 

The use of the mapping (6) can lead to a dramatic reduction in the error when 
there are multiple domains. In this case the overall accuracy of the Chebyshev 
interpolation is not sensitive to the location of the interface as long as the domains 
are such that there are not widely separated regions of rapid variation within each 
domain. In the solution of partial differential equations it may be important to 
determine the domains adaptively so that the regions of rapid variation are localized 
within each subdomain. The strategy of minimizing the maximum of (3) appears to 
accomplish this. Although we illustrate this behavior for only one example we have 
observed this behavior for other cases, not reported here, as well. 

The use of (7) leads to an improvement in accuracy over the unmapped Cheby- 
shev approximation but the improvement is not as great as in the case of the single 
domain. We note that the results obtained using (7) are slightly worse than in the 
unmapped case. This is due to the fact that the collocation points are nearly uni- 
formly spaced and since Q ^ 0, the spacing of the collocation points is greater in one 
of the subdomains. Setting Q = 0 we get errors very similar to the single domain 
case. 

In summary, results indicate that the use of the mapping (6) in conjunction with 
a multiple domain procedure can lead to a improvement in Chebyshev interpolation 
and pseudo-spectral methods for the case that the function to be approximated 
has multiple regions of rapid variation. Furthermore, this improvement can be 
obtained from adaptive procedures without prior knowledge of the location of the 
rapid variation. 

The mapping (6) can be extended to higher dimensions to deal with highly curved 
regions of rapid variation. One way this can be done is by allowing a to depend 
on the transverse coordinates similar to the procedure employed in [4]. Domain 
decomposition in several dimensions is usually based on rectangular domains. In 
this case such regions will not always be contained in one subdomain. 
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FIGURE la. (9) PLOTTED AGAINST THE INDEPENDENT VARIABLE X. 



FIGURE 2a. (10) PLOTTED AGAINST THE INDEPENDENT VARIABLE X. 




FIGURE 1b. (9) PLOTTED AGAINST THE TRANSFORMED VARIABLE S. FIGURE 2b. (10) PLOTTED AGAINST THE TRANSFORMED VARIABLE S. 
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